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THE UPPER ATMOSPHERE AS A 
MULTIPLE REFRACTIVE MEDIUM 
FOR NEUTRAL AIR MOTIONS 


H, Volland 


ABSTRACT 


Under the influence of gravity and heat conduction four plane characteristic 
waves obliquely incident on a horizontally stratified atmosphere can propagate, 
two of them upward and the two other downward. The two pairs of character- 
istic waves are the well known acoustic -gravity waves and the heat conduction 
waves. Molecular viscosity generates two further pairs of characteristic waves, 
the ordinary and the extraordinary viscosity waves. Ion drag and Coriolis force 
make the atmosphere anisotropic with respect to the characteristic waves. 
Their propagation characteristics for east to west and north to south propaga- 
tion differ from each other. 

Some analytical solutions of the eigenvalues of these eight characteristic 
waves are given in this paper. Numerical calculations of the eigenvalues de- 
pending on the parameters of the thermosphere, on frequency, azimuth and on 
angle of incidence are presented and discussed in some detail. 
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THE UPPER ATMOSPHERE AS A 
MULTIPLE REFRACTIVE MEDIUM 
FOR NEUTRAL AIR MOTIONS 


I. INTRODUCTION 

Since the work of Harris and Priester (1962) and of Hines (1960) it is well 
established that neutral air waves propagate within the thermosphere driven 
either by the solar EUV heat input within the thermosphere or originating in the 
lower atmosphere and propagating upward into the thermosphere. The paper of 
Harris and Priester stressed the importance of heat conductivity for the propa- 
gation characteristics of the diurnal waves. This is due to the fact that the ratio 


k/ p 


( k = coefficient of heat conductivity; p - mean density) is increasing with height 
because k depends only on temperature while p decreases exponentially with 
altitude . The change in the propagation behavior of gravity waves due to heat 
conduction and viscosity has been shown by Pitteway and Hines (1963). They 
calculated the eigenvalues of plane harmonic gravity waves and discussed the 
energy dissipation of these waves due to heat conduction and viscosity. Midgley 
and Liemohn (1966) made a full wave calculation of neutral atmospheric waves 
below 160 km altitude taking into account heat conduction and viscosity. They 
used a numerical technique to suppress the evanescent heat conduction waves 
and viscosity waves and calculated reflection coefficients of the lower atmos- 
phere with respect to gravity waves. 

More recent observations of traveling ionospheric disturbances (TID) made 
by Thome (1964) and Georges (1967) and observations of wave like structure in 
the thermospheric neutral density by Newton et al. (1968) show that gravity waves 
can reach altitudes as high as 500 km. At these heights however heat conduction 
waves and viscosity waves are no longer evanescent waves compared with gravity 
waves and may play a more or less important role. 

In order to deal with either the ray treatment or the full wave treatment of 
atmospheric waves one has to know the eigenvalues of the different wave modes. 
In this paper analytical and numerical calculations of the dispersion equation 
valid at thermospheric heights will be presented. For convenience we restrict 
ourselves to the simplest type of waves, namely plane harmonic nonducted 
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neutral atmospheric waves propagating obliquely through a horizontally stratified 
quiet atmosphere under the influence of heat conduction, viscosity, ion drag and 
Coriolis force. Under those conditions four pairs of characteristic waves exist: 

a. Acoustic-gravity waves 

b. Heat conduction waves 

c. Ordinary viscosity waves 

d. Extraordinary viscosity waves. 

Characteristic waves are only well defined within a homogeneous medium. 
The atmosphere however behaves like an inhomogeneous medium with respect to 
neutral air waves because the eigenvalues of the characteristic waves are func- 
tions of 


P P 

k = coefficient of heat conductivity 
rj = coefficient of viscosity 
p = mean pressure 

This leads to the construction of a model atmosphere consisting of a number 
of homogeneous slabs of constant temperature and constant parameters */p and 
• 77 /p. In each of these slabs the eight characteristic waves propagate uncoupled 
from each other. Coupling then occurs only at the boundary between adjacent 
slabs with different parameters. 

In this paper we shall solve the eigenvalue equation of the characteristic 
waves analytically and numerically within one homogeneous slab. We shall dis- 
cuss the behavior of the eigenvalues of the different wave modes as functions of 
angle of incidence, frequency and thermospheric parameters. The validity of 
such model atmosphere in full wave and ray treatments has been demonstrated 
in an additional paper (Volland, 1968). 
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II. THE BASIC EQUATIONS 

We start from the equations of conservation of mass, momentum and energy 
and the ideal gas equation of the neutral gas which are 

9/0 _ 

Xt + div (pv) - 0 

dv _ _ _ _ 

+ div cr + V p (v - Vj ) - 2Q x v p + grad p - pg 0 

( 2 ) 

dT . _ _ 

c v pyjT + P div v + (3 ~ div (k grad T) + u p v * (v - v i ) - 0 

p - /x p T = 0. 

p density; v = (u, v, w) velocity 
p pressure; T temperature 

v collision frequency between ions and one molecule 
cr viscous stress tensor 

v i 3 v k 2 1 

‘Xk - - + yy-~ 3 8 ik div v > element of cr 

t) coefficient of molecular viscosity 
0 earth's rotational vector 
g gravitational acceleration force 
c v specific heat at constant volume 

(v * B 0 ) B 0 

V- 

*o 2 


3 



ion velocity parallel to the geomagnetic field B 0 

3 Vj 
k 

i k 

viscosity heating 

« coefficient of heat conductivity 

/X = R/M 

R gas constant 
M molecular weight 

For solving the system of Equations (2) a consequent perturbation method is 
applied assuming that the time independent mean values like density p, pressure 
p and temperature T are already known and that the mean velocity of the air v 
is zero. The perturbation is considered to be a plane harmonic nondueted wave 
of angular frequency u> and of wave number k obliquely incident on a horizontally 
stratified plane atmosphere in which all parameters depend only on the vertical 
component z. Then all variables are functions of the coordinates x , y, z, and 
time t according to 

f(z) e ia>t-jk sin 0(x cos ^ + y sin ^ ) (2) 

6 is the angle of incidence and A is the azimuth of the plane of incidence of the 
wave with respect to geographic south. We assume that no internal energy 
source exists. 

The vertical component of energy transported by the different wave modes 
is continuous at any internal boundary. This implies the equivalence to Snell's 
law: 



k sin 6 - const. 


A - const. 


(4) 
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We normalize wavenumber k by an arbitrary constant wavenumber 



where 

C o = /tmT 0 (6) 


is the acoustic phase velocity at an arbitrary height z 0 . Since we exclude ducted 
wave modes, our model atmosphere is extended infinitely in vertical direction. 
Thus the normalized horizontal wavenumbers 



S 0 cos A o 



S 0 sin A 0 


can have all possible real values and remain constant throughout the whole at- 
mosphere. 

For convenience we rotate the coordinate system in such a manner that the 
new x-axis points in the A Q direction. Then the horizontal components of the 
velocity with respect to the new coordinate system are 


u - u cos A Q + v sin A 0 


v - - u sin A Q + v cos A 0 


( 7 ) 


The values Au, Av, Aw, Ap, Ap and AT are the deviations from the mean values 
and are considered to be small compared with C 0 , p (z), p (z) and T(z) of the 
quiet atmosphere. Thus all products of these small values can be neglected. 

We normalize these variables according to 
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Aw _ AT _ Au 

6 1 = C Q ; ^ = T 1 6 S = “q 1 ^ 


Ap _ «-AT' Au' 

“ r; e. - — — ; e, = : e 


'V/ 

Av 




( 8 ) 


The derivative with respect to z has been replaced by a prime (3/9 z = 1 ). 


Introducing the expressions of Equation (8) into the system of Equations (2) 
and considering the time and spatial dependence of Equation (3) the variables Aw" 
and A p can be eliminated which leads to a system of first order differential 
equations which in concise matrix form is 


I 


where 


e' “ jk 0 Ke - 0, 



(9) 


is the column matrix containing the 8 independent variables defined in Equa- 
tion (8). The coefficient matrix 



( 10 ) 


has the submatrices 

I 


6 



-1 


1 


0 




2 j (Aj - A 2 ) S - 2 j (A - A 2 S) 
0 2j (Aj - A) 
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0 


0 


0 


Rdi(l- jB 7 )- jS- 


jRBjd, 0 Rd^l- jB g )- jS 0 2 


Here the following abbreviations have been used: 


C Q - V7mT/ z = ! 


V + k 2 

x y 


C - UynT(z) 


2k 0 H 


2k 0 Hi 


2k 0 H 2 


7 

(y ~ i ) 


8 


s 


1 


v 


_3 j_R’ Z 1 co ’ Z 2 

4 y 


D 


1 


+ 


3 RZ 2 

4 7 


cos 2 




2 a; \ 

4Af*4A,A 2 -— j 



Zj sin $ sin A Q ± Z 2 sin I cos I cos A Q 



Z 1 sin t? cos A 0 ± Z 2 sin I cos I sin A Q 


Z x cos t) ± Z 2 cos 2 I cos A Q sin A Q 



cos I 


cos 


2 A, 


1 ~ cos 2 I 


V 


in 2 A„ 


sin 


in d 


sin 


+ 3 cos 2 d 


; sin I - 


2 cos d 




d geographical co -latitude 
I geomagnetic dipangle 
A Q azimuth 

In this calculation the geomagnetic field has been approximated by a dipole 
field with its axis parallel to the earth's rotational axis. Moreover, the derivative 


9 


of jx with respect to z has been neglected. But we have allowed for an altitude 
dependence of p, p, p and k . 

For comparison with earlier works (e.g., Hines, 1960; Eckart, 1960) we 
notice the following characteristic frequencies: 




(Brunt-Vaisala (11) 

Frequency) 


w h 


yg _ 

2 V 


dj G a 


Frequency co g is the critical frequency of gravity waves in a loss free iso- 
thermal atmosphere, cj is a stability parameter for small scale dynamics. 
o) h is a measure of the dissipation of wave energy due to heat conduction. The 
parameter V in w h has the dimension of a velocity and can be interpreted as the 
velocity of vertical heat transport within a static atmosphere in which the ratio 
«/ p remains constant (Volland, 1967). 

The parameters R and G are proportional to each other: 

2 f 

R = G (12) 
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with 


K 

f = ~ 2.5 

Vc v 


(Chapman and Cowling, 1959), indicating the inherent relationship between heat 
conduction and viscosity. R can be interpreted as some kind of Reynolds number 
because for acoustic waves the ratio between acceleration force and viscosity 
force in the equation of momentum is 
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V \ 2 

If we take into account a mean horizontal velocity of the air with components 
(u, v) then we only have to replace angular frequency co by an effective frequency 

W e f f ~ w ~ ^0 ^0 u 

where u is the component of mean velocity in the A 0 -direction of wave propaga- 
tion as can be derived from Equation (7). 


3 v 


n A v 


III. THE CHARACTERISTIC WAVES 

In order to solve Equation (9) uniquely it is necessary to know the conditions 
at the lower or the upper boundary of the model atmosphere. Physically appro- 
priate solutions require separation between waves transporting energy upward 
and waves transporting energy downward. Such waves are the characteristic 
waves or wave modes. We find them by transforming the matrix of the physical 
parameters e in Equation (9) 


e = P c (13) 

where the column matrix 



contains four upgoing waves a. and four downgoing waves b^ Introducing 
Equation (13) into Equation (9) gives 

c' = jk Q Nc - P _ 1 P' c (15) 

P in Equation (13) has to be chosen in such a manner that the matrix 
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in Equation (15) is a diagonal matrix. The elements in N are the eigenvalues 
of K defined by the eigenvalue equation 

| K + E | = 0 (17) 

(E unit matrix). 

Within a homogeneous medium in which the elements of K and P are con- 
stant the characteristic waves now have the solutions 



_ " i ko (z~ *o) 

‘“Oi e 


(18) 


and are independent of each other. However, the parameters G, R and Z 2 
in K are altitude dependent even within an isothermal atmosphere. The atmos- 
phere therefore behaves like an inhomogeneous medium with respect to neutral 
air waves. The matrices K and P are functions of height. The matrix P -1 P’ 
in Equation (15) therefore couples the different characteristic waves with each 
other at any height. 

If this coupling is weak, such that the elements of P _1 P’ are small com- 
pared with an approximate solution of Equation (15) is 


12 


or 


a. 

i 


b. 


r e - j ko j; 0 A i ( ^ )d e 

C 0i 


(19) 


This is the ray solution of Equation (9) (e.g., Budden, 1961). It means that the 
characteristic waves propagate uncoupled from each other through the inhomo- 
geneous atmosphere and that their amplitude and phase behavior is governed 
only by their eigenvalues k in the phase integral . 

The WKB solution as a second order approximation is defined by a special 
normalization of the elements of P such that the matrix P _1 P' has only zero 
elements in its diagonal. As a consequence the amplitude factors c 0i in Equa- 
tion (19) become height dependent. 

It can be shown by a numerical full wave treatment of Equation (9) that 
Equation (19) is in fact an excellent approximation at least for ascending gravity 
waves in thermospheric heights (Volland, 1968). 


IV. THE ISOTHERMAL MODEL ATMOSPHERE 

In view of the important role of the eigenvalues in any wave treatment our 
main task in this paper is to find analytical and numerical solutions of the eigen- 
value Equation (17). Special analytical solutions have the advantage that they 
allow a qualitative discussion of the behavior of the characteristic waves and 
that they help to identify the various wave modes in numerical treatments of the 
eigenvalue problem. 

For convenience we confine ourselves to an isothermal atmosphere. This 
has the advantage that it substantially simplifies the eigenvalue Equation (17). 
Moreover it leads to the classical formula of Hines (1960) for an adiabatic atmos- 
phere (k = rj = v = 0). This assumption is not really a restriction because in 
numerical full wave treatments one generally approximates the realistic atmos- 
phere by a number of isothermal slabs (e.g., Midgley and Liehmon, 1966; Volland, 
1968). A thickness of 1 km per slab is in most cases a sufficient approximation. 
In a ray treatment like Equation (19) a change in temperature can be taken into 
account by a height dependent parameter A(z). 
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A characteristic wave defined in the foregoing section is completely 
uncoupled from the other characteristic waves within one slab if the parameters 
of the coefficient matrix K are constant [see Equation (18)]. Then coupling 
takes place only at the boundaries to the adjacent slabs which have other parame- 
ters. If the slab is isothermal with exponentially decreasing pressure and 
density it follows from the definitions of G, R and Z 2 that the ratios 

— JL A — 

— i and 

p p ^ 

must be kept constant within the slab. This involves a constant collision fre- 
quency v and an exponential decrease of k and 77 with height such that 

A 2 = Aj = A = const. 

2 f 

R - — G - const. 

y 

d ! = I- 


Our model atmosphere therefore has an altitude dependence of the physical 
parameters as shown in Figure 1. In order to meet the dynamic boundary con- 
dition in hydrodynamics the pressure p(z n ) must be continuous _at the boundary 
between two adjacent slabs with different temperatures T n and T n + 1 . 

A more realistic model of the atmosphere could be found by approximating 
the real temperature profile within one slab by an exponential law 


T* 

n 

V - = const . 


This would lead to continuous pressure, density and temperature values at the 
boundaries and an k -dependence like p. It would greatly increase however the 
difficulties in solving Equation (17) and gives no real improvement in numerical 
calculations. 

The profiles shown in Figure 1 approach realistic atmospheric profiles if 
the thickness Az n of the slabs is sufficiently small. Numerical full wave and 
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ray optics calculations completely justify this kind of approximation including 
even the strange zigzag profiles of the coefficients of heat conduction « and 
viscosity 77 (Volland, 1968). 

V. EIGENVALUES OF AN ISOTHERMAL ATMOSPHERE 

From the theory of gravity waves within an isothermal nondissipative at- 
mosphere it is well known (Hines, 1960) that an upgoing wave increases its 
amplitude like 


e z/2H 


( 20 ) 


in order to keep the wave energy constant in an atmosphere where the pressure 
decreases exponentially with height. In a dissipative atmosphere the amplitude 
of a wave cannot surpass the amplitude of the equivalent nondissipative wave 
because part of the wave energy is transferred into internal energy of the sur- 
rounding gas by heat conduction, viscosity or collisions with the plasma com- 
ponent.* We therefore split from the eigenvalue a positive imaginary part j A: 

= q. + jA (21) 


which describes the amplification of the amplitude as shown in Equation (20). 

The remaining term q lt which we shall refer to as the eigenvalue, is com- 
plex: 


% 


a i ‘ iA 


( 22 ) 


Here the real part a ; is a measure of phase velocity V p : 




(n i = real refractive index) 


( 23 ) 


* 


This transfer of energy is a nonlinear process which is not included in our linearized theory. 
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and of the vertical component of group velocity 


V 

B 1 


d ( coa i ) 


a. 


1 a> 


- - S 


0 3S„ 


( 24 ) 


The angle between the vertical and the phase velocity is 



arc tg 


a. 


(25) 


while the angle between the vertical and the direction of the group velocity within 
the (a 4 , S 0 ) — plane is defined as 


| = - arc t g ^ - (26) 

The imaginary part - j /L is responsible for the dissipation of wave energy. We 
define as upgoing or downgoing waves 


/3j ^ 0 (upgoing) 

£ 0 (downgoing) 


(27) 


where the equal signs only stand for nondissipative waves {« = rj - v - 0). 


VI. CLASSIFICATION OF THE CHARACTERISTIC WAVES 
IN AN ISOTHERMAL ATMOSPHERE 

In this section we shall solve the eigenvalue Equation (17) analytically 
for some very special cases in an isothermal atmosphere. We shall define 
the different characteristic waves — or wave modes. As follows immedi- 
ately from the kind and number of the system of differential Equations (9) 
there exist four pairs of wave modes, each pair consisting of an upgoing 
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and a downgoing wave. These four wave pairs are 


Acoustic-gravity waves (q r q 2 ) 
Heat conduction waves (q 3 , q 4 ) 
Ordinary viscosity waves (q 5 , qg) 


Extraordinary viscosity waves (q 7 , q 8 ) 

In general all eight eigenvalues q. are different from each other showing 
that waves belonging to one pair may have different propagation characteristics. 
In our special analytical solutions we shall find pairs of eigenvalues which only 
differ in sign (e.g., q x = -q 2 ,etc.). 

a. Extraordinary Viscosity Waves 


We start with the extraordinary viscosity mode because this mode can be 
easily separated from the other modes. From the eigenvalue Equation (17) it 
can be seen that this wave mode is coupled with the other waves only through 
Coriolis force and (or) ion collisions. If Z x , Z 2 are small compared with 1 
the parameters B n in K are small and can be neglected. Then the eigenvalue 
Equation (17) has the form 
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and 


K. + X. E = 

6 1 


^ R - i s « 



Now the eigenvalue equation can be separated into two determinants: 


I K s + E | • | K g + X. E | = 0 

which implies: 

I K g + E | - 0 

| K 6 + X. E I =0 


( 30 ) 


(31) 

(32) 


The solution of Equation (32) is 


\ 7 - j A - q 7 - + j|^A 2 + S Q 2 + j R (33) 

8 8 

(Zj = Z 2 = 0) 

Q 7 are the eigenvalues of the pair of extraordinary viscosity waves. The upper 
sign in Equation (33) stands for the ascending wave, the lower sign stands for the 
descending wave as one can check immediately from the definitions given in 
Equations (27). This wave mode is a purely transverse wave inasmuch as only 
the horizontal component Av orthogonal to the direction of wave propagation is 
involved. It is heavily attenuated for large values of R and becomes a purely 
evanescent wave at R = 0. 

Equation (33) is also exact for east-west propagation (A 0 = ±90°) at the 
equator ( 0 - 90°) even for arbitrary values of Z x and Z 2 . In this case the wave 
is likewise completely uncoupled from the other waves. 
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From Equations (23) and (24) we derive phase and group velocity at vertical 
incidence as 


r 

± 


2C 0 A 

R 


r 


A >> R 



for -l 


A << R 


r 


\ - 1 


2 V 

p 




A >> R 

for < 


A << R 




(34) 


Phase and group velocity have the same direction, 
b. Ordinary Viscosity Waves 

The remaining Equation (31) can similarly be separated for vertical inci- 
dence S 0 = 0. Then Equation (31) becomes 



+ Ml 


K i + M 


o 


o 


= 0 


K t + M 


(35) 


K ? looks exactly like K 6 in Equation (32) for S Q = 0 and gives the eigenvalues of 
the pair of ordinary viscosity waves 


q s - + j l/A 2 + j R 

6 


( 36 ) 


(S 0 = = Z 2 = 0) 
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which have the same behavior as the eigenvalues of the extraordinary viscosity 
waves. 


c. Acoustic -Gravity Waves and Heat Conduction Waves 


The eigenvalues of the remaining matrix K x in Equation (35) can be found 


for S 0 = 0 as 


2 


>= +i 


% 

4 


jj/vTT^T^ ?i|( 


(G + R*) 


_ 


4GR' 


1 - 


2 j G 

yd 


(37) 


(S 0 = Zt = Z 2 = 0) 


with 


r 


R* = -y (1 - S) 




yR 


R << 1 




(y - i)y 

2 f 


for 


R R >> 1 


The negative sign within the square root in Equation (37) is related to acoustic- 
gravity waves, the positive sign is related to heat conduction waves. The signs 
outside the square root are again due to ascending and descending waves. 


Approximate solutions of Equation (37) are 


*\ 


r 


4 * *Ti 




for G, R << 1 


+ jVA 2 + 


'J 


2 ) G 

y 
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±/l - A 2 


2 




r 




for G, R >> 1 


q 3 

4 


J 


Tj^A 2 + 2 j G 


( 38 ) 


(S 0 = Zj = z 2 = 0) 


Equations (38) indicate as expected that for small Reynolds numbers (R <" < 1 ) 
all wave modes including the acoustic-gravity mode are evanescent waves be- 
having like viscosity waves while for large Reynolds numbers (R > > 1) heat 
conduction waves behave like viscosity waves but acoustic -gravity waves are 
propagation waves for A < 1 (acoustic wave range) and are evanescent waves 
(at vertical incidence) for A > 1 (gravity wave range). 


VII. APPROXIMATE SOLUTIONS OF THE EIGENVALUE EQUATION 
AT LARGE REYNOLDS NUMBERS (77 - 0) 

An important approximation of the general eigenvalue Equation (17) can be 
found if one makes the assumption 


R 



(39) 


In this case we again can split the determinant |K + X. E| into two nearly inde- 
pendent determinants: 


+ x iE 


K + X. E 

1 


0 (40) 


K 4 + \ E 
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We do this by multiplying the 5th and the 7th column of |k + Xj E| with factors 
h x and h 2 , respectively, and add both columns to the first column. We next 
multiply the 5th and 7th column by 

s o S„ 

~~~ h, and h. 


and add both columns to the second column. Here it is 


hj (B 3 Bg + B 2 B g + j B 2 )/A 

h 2 = (B 3 B 7 - B 2 B 5 + jB 3 )/A 

h 3 - ( 1 “ 3 B g )/A (4D 

h 4 = jB s /A 

A = 1 - B 5 B 6 - B 7 B 8 - j(B 7 + B 8 ) 

This procedure eliminates the R dependence of the off diagonal submatrices in 
Equation (17), since S 0 and 


4 


3y 


Their elements therefore become negligibly small compared with the elements 
of the submatrix K 4 + X . E if R goes to infinity. Thus the eigenvalues of the 
viscosity waves can be determined from the determinant 

|K 4 + X. E| = 0 


as 
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+ jy A 2 + s 0 2 + *9 (B ? + Bg) + j R ± jRVB 5 B 6 - "T(B 7 - B 8 )' 


The remaining determinant in Equation (40) has the form 


^ + hj S 0 - 2 j A - 1 ; 


B 3 S o 


|Kj +X..EI 


^i-h 6 S 0 


d 2 l 1 + 2jG 


h s ~ i ( Z 2 COs2 1 ~ B l h l “ B 4 B 2 ) 


h 6 " i( B ih 3 + B 4 h 4 ) 


It can be shown that the element 


4 (y-l) 

K = — -v — ~ 0 4 

"•24 3 ' f 
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which is the only direct contribution of R via the parameter S in this determinant 
can entirely be neglected as long as 

A » K 24 . 


This is in the gravity wave range. 

The negligence of the term K 24 in Equation (43) means that we treat the 
system of Equations (2) without taking into account viscosity (17 - 0). The 
transformation of the determinant Equation (40) has then the effect of eliminating 
the horizontal components of the velocity in the momentum Equations (2) which 
are coupled via Coriolis force and (or) ion drag with the vertical component of 
the velocity Aw and with the pressure Ap. The terms hj and h 6 in Equation (43) 
are now contributions of the vertical wind velocity and of the horizontal pressure 
gradient. They are of the order £ 1 and therefore small compared with A or 

!M- 

If we neglect in Equation (43) the element K 24 and all terms which include 
the difference hj - h 6 , the solution of Equation (43) becomes 


a 


<*1 

2 




+ i ^ 2 + s o 2 + j G + 2j 


G ( h 8 + s 0 2B 2 h 3 ; 


^3 


j 


(v - 0; hj - h 6 - 0) 


( 44 ) 


with 

h ? = Sp (1 + hjh 6 ) + (l-h 5 ) (y-S 0 2 h 3 >- jS 0 A(h 1 + h 6 ) 

2jS 0 A 

h g - h 7 - (y- 1) (1 - h s ) + — (h 1+ h 6 ) 

2 

d 2 y co 2 
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Again the upper sign within the square root is for acoustic -gravity waves, the 
lower sign is for heat conduction waves. 

For east- west propagation (A Q = ±90°) at the equator (d = 90°) it is hj = h 6 , 
and Equation (44) is exact provided K 24 = 0. 

For Zj = Z 2 = 0 it is 



Then Equation (44) becomes identical with Pitteway and Hines (1963) Equation (36). 

From our assumption made in Equation (39) we should expect that the ap- 
proximate Equation (44) breaks down for small values of R (and therefore of G). 
Surprisingly enough we shall find from numerical calculations given in the next 
section that Equation (44) is a reasonable approximation even for 


G > A > 1 


which is within the range of gravity waves as long as S 0 is not too large. 

VIII. NUMERICAL SOLUTIONS OF THE EIGENVALUE EQUATION 

Apart from some special analytical solutions which are given in the previous 
sections, Equation (17) has to be treated by numerical methods. A convenient 
program available in the SHARE Program Catalog is SAD 3099 (EIG 4) which 
solves eigenvalues of complex matrices and is coded for the IBM 7090/94. This 
program has been used for the following calculations . 

The atmospheric model used is the Harris -Priester model 5 at 1200 local 
time in 200 and 400 km height (CIRA, 1965). The parameters of this model to- 
gether with additional numerical values are given in Table I. 
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Table I 


Numerical Data Used in the Calculations of Eigenvalues 


z (Km) 

"aO" 1 ) 

a, h( s_1 ) 

V ( S l ) 

fi(s 1 ) 

7 

f 

! 

C o ( m /s) 

V (m/s ) 

200 

8.93 x 10" 3 

1.09 

1. x 10“ 4 

7.27 x 10~ 5 

1.5 

2.5 

773.0 

6.35 

400 

6.51 x 10~ 3 

3.92 x 10' 2 

1. X 10" 4 

7.27 x 10 -5 

1.5 

2.5 

999.0 

166.0 


In each of the following figures the eigenvalues of the different wave modes 
has been plotted versus the normalized horizontal wave number 


k 2 

y 




The curve parameter is the angular frequency co (in sec" 1 ). The left side of each 
figure gives the real part a i of the eigenvalue. The right side gives its negative 
imaginary part multiplied by k 0 : k^. We restrict ourselves to ascending 
waves. The eigenvalues of descending waves are not very much different from 
the eigenvalues of the equivalent ascending waves (apart of course from the sign). 

Figures 2-5 show results of numerical claculations of the eigenvalue Equa- 
tion (17) (full lines). For comparison the dashed lines have been calculated from 
the approximate formula Equation (44) with the same numerical data. The propa- 
gation conditions in these calculations are east to west propagation (A 0 - 90°) at 
the equator (tf = 90°). In order to arrange that the various curves fit within one 
figure, the ordinate of ctj in Figures 2 and 6 has been divided into two different 
scales — a linear one for a t > - 1, and a logarithmic one for a t < - 1. 

From the definitions of magnitude and direction of phase velocity [Equations 
(23) and (25)] it follows that the vector from the origin to a point of the dispersion 
curves gives the ratio between velocity of sound to phase velocity and the direc- 
tion of phase velocity. We observe in Figure 2 the well known downward directed 
phase velocity of the ascending gravity waves {co < co a ) while their energy propa- 
gation vector defined by Equations (24) and (26) is orthogonal to the dispersion 
curves of a. (S 0 ) and is directed upward. 
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HEAT CONDUCTION WAVE 
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ORDINARY VISCOSITY WAVES 



K 


Figure 4. Real part a 5 (left) and negative imaginary part k Q /3j (right) of eigenvalues of ascending ordinary viscosity 

For details see text of Figure 2. 









ACOUSTIC-GRAVITY WAVES 
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Figure 6. Real part a.] (left) and negative imaginary part k Q 5-j (right) of eigenvalues of ascending acoustic-gravity 

in 400 km height. For further details see text of Figure 2. 




Acoustic waves (co > co a ) with Reynolds numbers R £ 1 show a normal 
behavior, but become evanescent waves at R < 10" 1 as expected from Equa- 
tion (38). 

The intermittent frequency (o> a ~ a = 10" 2 sec" 1 ) shows acoustic wave 
features at S 0 < 1 and gravity wave features at S 0 > 1 with a continuous transi- 
tion region. From Figure 2 (attenuation factors k 0 /3 x ) it follows however that 
the gravity wave range of this frequency is heavily attenuated. 

From the definition of the attenuation factor [Equations (22) and (27)] we 
can deduce that the wave amplitude of an upgoing wave remains constant if 


A 


1 

2H 


This is the general behavior for gravity waves (co < coj and S 0 < 10 in Figure 2. 
The frequency 10" 3 s 1 has a minimum of attenuation at S Q = 1.5. At this 
frequency and horizontal wave number (k ~ 2 x 10" 3 km" ! ) we expect therefore 
especially good propagation conditions in this altitude range. The atmosphere 
seems to behave like a frequency and height dependent selective filter with re- 
spect to gravity waves. 

Acoustic waves with frequencies co > Is -1 are heavily damped and do not 
contribute to energy transport in this height. 

The dashed lines in Figure 2 calculated from the approximate formula 
Equation (44) indicate that this formula is a reasonable approximation within the 
gravity wave range (A > 1 or co a > co). It breaks down in the acoustic range. 
There one expects from Equation (44) a phase velocity 


V 

p 



7 


for co -> oo while Equation (38) predicts V p -» 0. 

Figures 3 to 5 show real and imaginary part of eigenvalues of heat conduc- 
tion wave, ordinary and extraordinary viscosity waves. All three wave modes 
behave very similar in the whole frequency range. Heat conduction waves are 
slightly less attenuated than viscosity waves 
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For A > 1 and S Q > 10 real and imaginary part of the eigenvalues in 
Figures 3 to 5 are of the same order of magnitude as the gravity waves (Figure 2). 
The sign of a i in Figures 3 to 5 is always positive at this height indicating that 
energy and phase normal are equally directed. 

Within the acoustic range (A < 1) all three wave modes in Figures 3 to 5 tend 
to become heavily attenuated evanescent modes. ( k 0 /3. = 1 is equivalent to a 
decrease in amplitude of 1/e after a vertical propagation path of 1 km.) Neither 
mode ever shows an increase in amplitude with height. 

The dashed lines in Figures 3 to 5 have been calculated from Equations (42) 
and (44) and give evidence of the rather good approximation which these formulae 
provide . 

In order to check the validity of Equation (44) as compared with the exact 
Equation (17) the calculations have been repeated in Figure 6 for gravity waves 
using the parameters at a height of 400 km (see Table I) and the same propagation 
conditions (A 0 = 90°; d = 90°). We observe again a reasonably good agreement 
between the exact solutions from Equation (17) (full lines) and the approximative 
solution from Equation (44) (dashed lines) for the gravity wave range (A > 1). 

This agreement breaks down for A > G or for very large S 0 values. A ~ G 
holds in about 500 km. Below this height the error arising from Equation (44) 
seldom exceeds 50% in the gravity wave range. In view of the great simplifica- 
tion given by Equation (44) as compared with the exact expression Equation (17) 
this error seems to be of an allowed order. However the approximate decoupling 
of the viscosity waves from the gravity and heat conduction waves greatly sim- 
plifies any full wave treatment. It justifies for instance Harris' and Priester's 
(1962) calculations in which the coefficient of viscosity has been neglected. 

In order to investigate the influence of Coriolis force and ion drag we com- 
pare in Figures 7 to 10 the real and imaginary part of eigenvalues of gravity 
waves at 200 km for south to north propagation ( A 0 = 180°) with east to west 
propagation (A Q = 90°) (dashed lines) at the equator (d = 90°) (Figures 7 and 8) 
and at 45° northern latitude (d = 45°) (Figures 9 and 10). Moreover we compare 
these dispersion curves with the eigenvalues of gravity waves where Coriolis 
force and ion drag have been neglected (Zj = Z 2 = 0) (full lines). Again the exact 
expression Equation (17) (Figures 7 and 9) and the approximate solution Equa- 
tion (44) (Figures 8 and 10) have been used. 

From these figures we see that Coriolis force and ion drag makes the at- 
mosphere anisotropical with respect to gravity waves. In Figure 7 the eigen- 
values for Zj = Z 2 =0 (full lines) are equal to the eigenvalues at south to north 
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propagation (A 0 = 180°). Coriolis force and ion drag therefore only affect east 
to west propagation (A = 90°) (dashed lines). The attenuation minimum at S 0 ~ 1 
becomes deeper for south to north paths as compared with east to west paths 
suggesting more favorable propagation conditions for south to north paths. 

In these and all the following calculations it has been found that for west to 
east propagation the eigenvalues do not differ more than 10% from the values of 
east to west propagation. The same is true for north to south paths in relation 
to south to north paths. 

Figure 8 gives the same results as in Figure 7 but calculated from Equa- 
tion (44) {rj = 0). For comparison the circles are drawn from the full lines in 
Figure 7 giving the exact calculation for Z t - Z 2 = 0. The dotted lines have 
been calculated from the conditions Z 2 = 0, A 0 = 90° indicating small influence 
of Coriolis force at the equator. From comparing Figures 7 and 8 we 
again observe a rather good agreement between values calculated from Equa- 
tions (17) and (44). We note that even special features of the single curves in 
Figure 7 are truly reflected in Figure 8. 


Finally Figures 9 and 10 give results at 45° northern latitude (d = 45°). 

Here we observe for the first time a remarkable difference between the results 
from Equations (17) and (44). It arises for the case co = 10 -s sec -1 , negligible 
ion drag (Z 2 = 0) and S 0 > 6. Here the exact a j -curves (dotted curves) bend 
toward zero and then become positive for increasing S 0 values. The equivalent 
a x curve in Figure 10 (dotted curve) behaves in the normal manner and shows 
no difference for A Q = 180° and A 0 = 90°. It differs however appreciably from 
the full curve (Z 1 = Z 2 = 0) giving evidence of the large influence of Coriolis 
force on the propagation characteristics in medium latitudes. If we add ion drag 
which in our calculations is of the same order of magnitude as the Coriolis force 
we find the dashed curves in Figures 9 and 10. Here ion drag has the tendency 
to reduce the influence of the Coriolis force. 

At co = 10" 4 sec" 1 the attenuation factor is larger in Figure 10 than in 
Figure 9 for Z 2 = 0. Again ion drag reduces substantially this difference. 


IX. CONCLUDING REMARKS 

It has been shown that eight obliquely incident plane characteristic waves 
can propagate through a horizontally stratified atmosphere, four of them 
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ascending and the other four descending. The four pairs of characteristic waves 
are the well known acoustic-gravity waves, the heat conduction waves and ordi- 
nary and extraordinary viscosity waves. Heat conduction waves and viscosity 
waves are named according to their relation to finite heat conductivity and finite 
molecular viscosity. 

Analytical solutions of the eigenvalue equation have been given which allow 
one to identify the different characteristic waves. 

Heal and imaginary parts of the eigenvalues of the characteristic wave 
modes have been calculated numerically. Through the use of these eigenvalues 
the general behavior of the different wave modes has been discussed in detail. 

The main results are the following: 

a. Under the influence of heat conduction, viscosity and ion drag acoustic- 
gravity waves become dissipative waves. Part of the wave energy is 
transferred into internal energy of the surrounding gas. This energy 
dissipation depends on angle of incidence, parameters of the atmosphere 
and frequency. In the gravity wave range there exists an attenuation 
minimum at a certain angle of incidence . The atmosphere therefore 
behaves like a frequency and height dependent selective filter with re- 
spect to gravity waves. 

b. Under the influence of Coriolis force and ion drag the atmosphere be- 
haves like an anisotropic medium. East to west propagation charac- 
teristics of gravity waves differ from north to south propagation char- 
acteristics. Ion drag tends to reduce the influence of the Coriolis force. 

c. Within the gravity wave range and at altitudes below 500 km viscosity 
waves are nearly uncoupled from gravity waves and heat conduction 
waves. The coefficient of viscosity therefore can be neglected in the 
treatment of gravity waves at thermospheric heights. 
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